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Quantification  of  macrobenthic  effects  on  diagenesis  using  a  multicomponent  inverse 
model  in  salt  marsh  sediments 

Yoko  Furukawa 1 

Naval  Research  Laboratory,  Seafloor  Sciences  Branch,  Stennis  Space  Center,  Mississippi  39529 


April  C.  Smith  and  Joel  E.  Kostka 

Department  of  Oceanography,  Florida  State  University,  Tallahassee,  Florida  32306 

Janet  Watkins 

Naval  Research  Laboratory,  Seafloor  Sciences  Branch,  Stennis  Space  Center,  Mississippi  39529 

Clark  R.  Alexander 

Skidaway  Institute  of  Oceanography,  Savannah,  Georgia  31411 
Abstract 

Using  a  multicomponent  inverse  model,  we  quantified  the  rates  of  organic  matter  (OM)  remineralization  and  the 
relative  importance  of  major  terminal  electron  acceptors  [Fe(III)-(oxy)hydroxides  and  SOJ  ]  in  salt  marsh  sediments 
with  varying  degrees  of  bioturbation  and  vegetation  at  Skidaway  Island,  Georgia.  The  model  determined  the  rates 
of  OM  diagenesis  by  seeking  simultaneous  agreement  between  measured  and  model-calculated  depth -concentration 
profiles  of  multiple  major  redox  species  while  using  the  biological  transport  parameters  determined  from  direct  field 
observations.  The  OM  degradation  rates  are  found  to  be  greater  and  penetrate  deeper  in  sediments  with  vegetation 
or  bioturbation  than  in  sediments  with  limited  macrobenthic  components,  in  which  organic  matter  degradation  is 
restricted  to  the  immediate  vicinity  of  water-sediment  interface.  The  results  also  confirmed  previous  observations 
that  Fe(III)-(oxy)hydroxides  are  much  more  important  than  sulfate  as  terminal  electron  acceptors  in  the  heavily 
bioturbated  station  and  provided  the  numerical  values  to  the  depth-dependent  relative  importance.  Solid-phase  Fe 
is  recycled  as  a  terminal  electron  acceptor  ^30  times  in  the  bioturbated  station,  because  biological  mixing  repeatedly 
moves  reduced  Fe(II)  back  into  the  aerobic  environment,  where  it  is  reoxidized  to  form  Fe(III)-(oxy)hydroxides 
and  is  reused  as  a  terminal  electron  acceptor.  Vegetation  appears  to  have  little  influence  on  net  solute  transport,  but 
it  stimulates  microbial  activities  and  significantly  enhances  the  remineralization  of  OM  at  depths. 


The  salt  marsh  is  a  complex  geochemical  system  with 
highly  coupled,  simultaneous  processes  such  as  primary  pro¬ 
duction,  microbial  degradation  of  organic  matter  (OM),  and 
macrobenthos  feeding.  Along  the  east  coast  of  the  United 
States,  salt  marsh  primary  production  is  dominated  by  a  spe¬ 
cies  of  smooth  cordgrass,  Spartina  alterniflora,  whose  root 
systems  are  the  conduits  for  chemical  exchange  between  sed¬ 
iments  and  overlying  water  or  air,  as  well  as  the  primary 
source  of  OM  to  sedimentary  microbial  communities  (Chil¬ 
ders  et  al.  1993).  Salt  marsh  sediments  are  often  heavily 
bioturbated  by  fiddler  crabs,  Uca  spp.  (Teal  1958).  These 
macrobenthos  are  directly  involved  in  chemical  mass  trans¬ 
fer  by  producing  and  metabolizing  OM  and  mobilizing 
chemicals  (Kristensen  et  al.  1991;  Nithart  et  al.  1999),  while 
indirectly  affecting  the  course  of  diagenesis  by  modifying 
microbial  habitats  (Goni-Urriza  et  al.  1999;  Marinelli  et  al. 
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2002).  The  complete  characterization  of  salt  marsh  system 
benefits  from  a  comprehensive  approach  in  which  all  critical 
processes  are  considered  concurrently  in  a  coupled  manner. 
Computational  reactive  transport  modeling  has  become  an 
effective  tool  for  the  comprehensive  descriptions  of  aquatic 
sedimentary  systems  undergoing  early  diagenesis. 

A  unique  utility  of  properly  constrained  modeling  is  the 
ability  to  quantify  properties  and  processes  that  cannot  be 
directly  measured,  such  as  the  in  situ  rates  of  OM  degra¬ 
dation  and  terminal  electron  acceptor  (TEA)  reduction.  In 
diqgenetic  modeling,  forward  modeling  proceeds  by  inter¬ 
actively  comparing  model-calculated  concentration  profiles 
of  geochemical  species  to  measured  profiles,  whereas  inverse 
modeling  uses  measured  profiles  as  a  direct  forcing  to  de¬ 
termine  hard-to-measure  properties  such  as  reaction  rates 
and  transport  coefficients.  Using  directly  measured  or  esti¬ 
mated  forcing  variables  (e.g.;  sediment  accumulation  rates, 
biological  reworking  depth,  and  overlying  water  chemistry), 
previous  forward  modeling  studies  have  hind-casted  the 
rates  of  Fe  and  Mn  cycling  in  shelf  sediments  (Van  Cappel- 
len  and  Wang  1996),  rates  and  mechanisms  of  denitrification 
at  the  global  scale  (Middelburg  et  al.  1996),  and  redox  os¬ 
cillation  along  infauna  burrows  (Furukawa  2001).  The  in¬ 
verse  approach  has  been  used  to  determine  consumption 
rates  of  major  redox  species  (Berg  et  al.  1998),  and  bioir¬ 
rigation  coefficients  (Meile  et  al.  2001). 
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Table  1.  Geochemical  species  and  reactions  considered  in  this  study 


Porewater  species 
02 
NO^ 
nh4+ 

SOf 

TS 

Fe2+ 

SC02 

Particulate  species 
OM 

Fe(OH)3 

FeS 

NH4(ads) 


Oxygen 

Nitrate 

Ammonium 

Sulfate 

Aqueous  sulfide  species  (=H2S  +  HS") 

Ferrous  ion 

H2C05  4-  HCOf  +  CO|- 

Organic  matter  with  the  model  formula  of  (CH20)(NH;)v/„  where  x:y  is  the  C:N  ratio 
Poorly  crystalline  Fe(III)-(oxy)hydroxides  expressed  as  Fe(OH)3 
Precipitated  ferrous  sulfides 
Surface-adsorbed  ammonium 


Primary  redox  reactions: 


OM  +  o2  ->  xco2  +  ^nh;  +  H20 

/?j  (mol  dm“3-solids  s~‘) 

OM  +  O.8NO3-  +  0.8H*  ->XC02  +  +  0.4N2  +  1.4H20 

R2  (mol  dm_3-solids  s~l) 

OM  +  4Fe(OH)3  +  8H+  ->  SCO,  +  ~NH4+  +  4Fe2+  +  11H.0 

x  2 

(mol  dm‘3-solids  s_1) 

OM  +  Q.5S0J-  +  H+  ->X0O2  +  -NH4+  +  11H20  +  0.5H,S 

X 

R4  (mol  drrr3~solids  s-') 

Secondary  reactions: 

Fe2+  +  H3,  +  ^H20  ->  Fe(OH),  +  2H+ 

Rates: 

Rs  (mol  L“ ‘-porewater  s-r 

NH4+  +  202  ->  N03-  +  HjO  +  2H+ 

H2S  +  202  -¥  SOr  +  2H+ 

“H,S  +  Fe(OH),  +  2H+  -4s°  +  Fe2+  +  3H20 
^  z 

R6  (mol  L“’-porewater  s”r 
R7  (mol  L~ ‘-porewater  s-r 

Rs  (mol  drrr3-so1ids  s_1) 

FeS  +  202  -»  Fe2t  +  SO2" 

Fe2+  +  H2S  ->  FeS  +  2H' 

R<j  (mol  dm~3-solids  s_1) 
Rw  (mol  L_,-porewater  s_1 

Equilibrium  reactions: 

NH|  <->NH4+(ads) 

Equilibrium  constant: 

kN 

Recent  studies  in  salt  marshes  of  the  U.S.  Atlantic  coast 
have  documented  the  significant  recycling  of  iron  and  sulfur, 
and  complex  spatial  overlapping  of  the  TEA  sequences  (Lu¬ 
ther  et  al.  1986;  Kostka  et  ai.  2002#;  Gribshoit  et  al.  2003; 
Koretsky  et  al.  2003).  This  recycling  has  been  attributed  to 
the  repeated  reoxidation  of  reduced  iron  and  sulfur  caused 
by  bioturbation,  bioirrigation,  and  phytoirrigation.  Because 
of  the  inherent  complexity  of  salt  marsh  systems,  it  is  dif¬ 
ficult  to  quantify  OM  remineralization  pathways  by  experi¬ 
mental  methods  alone.  A  modeling  approach  in  which  all 
critical  reactions  and  transport  processes  are  simultaneously 
considered  would  allow  comprehensive  and  quantitative  de¬ 
scriptions  of  the  complex  systems. 

The  objective  of  the  present  study  was  to  quantify  the 
rates  and  pathways  of  OM  remineralization  through  multi- 
component  inverse  modeling  constrained  by  extensive  field 
data.  The  realistic  model  descriptions  of  biological  transport 
processes  were  emphasized,  to  assess  their  effect  on  diage- 
netic  chemical  mass  transfer.  The  inverse  modeling  allowed 
the  extraction  of  hard-to-measure  parameters  using  directly 
measured  variables.  We  simultaneously  applied  inverse  mod¬ 
eling  to  concentration  profiles  of  multiple  species,  to  mini¬ 


mize  the  uncertainties  in  model  results  generated  by  the  in¬ 
dependent  uncertainties  in  rate  description  of  each  individual 
species. 

Materials  and  methods 

Reactions  and  rate  expressions — The  geochemical  species 
and  reactions  considered  in  the  study  are  defined  in  Table  1 . 
OM  degradation  pathways  in  the  upper  -0.5  m  of  salt  marsh 
sediments  have  been  considered  to  be  dominated  by  micro¬ 
bial  SO l~  respiration  (Howarth  and  Giblin  1983;  King 
1988),  whereas  recent  studies  have  shown  evidence  for  the 
importance  of  Fe(III)-(oxy)hydroxide  respiration  (Lowe  et 
al.  2000;  Kostka  et  al.  2002a).  Therefore,  we  investigated 
these,  in  addition  to  the  respiration  of  02  and  NOj,  as  the 
primary  redox  reactions.  The  reaction  expressions  are  based 
on  those  found  in  previous  work  (Froelich  et  al.  1979;  Van 
Cappellen  and  Wang  1996),  which  have  been  tested  in  the 
context  of  multicomponent  reactive  transport  models  (Boud¬ 
reau  19966;  Archer  et  al.  2002).  The  reactions  have  since 
been  updated  to  express  the  inorganic  nitrogen  product  of 
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aerobic  OM  degradation  to  be  NH4+  rather  than  N03  (Herbert 
1999;  Berg  et  al.  2003). 

On  the  basis  of  the  reactions  in  Table  1,  the  rates  of  pri¬ 
mary  and  secondary  reactions,  as  well  as  the  net  production 
rates  of  geochemical  species,  can  be  expressed  (Table  2). 
The  OM  degradation  rates  are  considered  to  follow  Monod 
kinetics  (Boudreau  1996a),  whereas  the  secondary  reactions 
are  presumed  to  follow  the  simple  bimolecular  rate  laws  that 
are  widely  used  in  multicomponent  reactive  transport  models 
(Van  Cappeilen  and  Wang  1996;  Berg  et  al.  2003). 

Mass  conservation  equations— The  multicomponent  in¬ 
verse  model  used  in  the  study  is  based  on  the  one-dimen¬ 
sional  mass  conservation  equation.  For  solute  species,  it  can 
be  written  as 


,  / wAjcY 

* =  co,hi^:j  ■ 


2D, 

wAx 


(4) 


(Boudreau  1996a).  C,  (mol  L~!)  and  D,  (m2  s’1)  denote  the 
concentration  and  net  diffusion  coefficient  (Ds  +  DB)  at  the 
zth  node,  and  Ax  is  the  distance  between  adjacent  nodes.  A 
similar  ODB  can  be  written  for  solid  species  by  partially 
discretizing  Eq.  2. 

By  assuming  that  the  geochemical  profiles  are  at  steady 
state,  the  left-hand  side  of  Eq.  3  is  set  at  dCJdt  =  0,  and 
the  right-hand  side  can  be  reorganized  as 

aaiCi-[  +  bbiC}  +  +  dd\  =  0  (5) 

where 


bc  =  ^ 

dt  (p  dx 


dC 

<p(Ds  +  DB)~ - cpxvC\ 

dx 


+  a(C0  -  C)  +  Xi? 

(1) 

or  under  the  no-compaction  (i.e.,  constant  porosity)  assump¬ 
tion, 


dC  =  8_ 
dt  dx 


C Ds  +  Dg)j- 
dx 


dC 

w—  +  a(C0  -  C)  +  ZR  (!') 


(Berner  1980;  Boudreau  1996a),  in  which  the  time  (/,  s)- 
dependent  concentration  of  a  dissolved  species  (C,  mol  L1 
pore  water)  along  the  vertical  axis  (x,  m)  is  determined  by 
its  diffusive  transport,  accumulation  due  to  sedimentation, 
transport  by  irrigation,  and  net  rate  of  production  and  con¬ 
sumption  reactions  &R,  mol  L'1  s'1).  As  (m2  s'1)  is  the 
tortuosity-corrected  molecular  diffusion  coefficient,  D„  (m2 
s_l)  is  the  biodiffusion  coefficient,  <p  is  the  porosity,  a  (m_1) 
is  the  bioirrigation  coefficient,  w  is  the  sediment  accumula¬ 
tion  rate  (m  s'1),  and  C0  (mol  L'1)  is  the  concentration  at 
water-sediment  interface  (WSI).  Similarly,  the  mass  conser¬ 
vation  equation  for  solid  species  can  be  written  as 

dC,  _  d  /  dC\  dC,  „ 

3/  "  +  (2) 

where  C,  is  the  concentration  of  the  solid  species  (mol  dm~3 
of  solids)  and  XRS  is  the  net  rate  of  production  (mol  dm’3 
solids  s’1).  Irrigation  and  molecular  diffusion  are  not  con¬ 
sidered  for  solids. 

The  partial  discretization  of  Eq.  1 '  for  the  zth  node  on  an 
even  grid  yields  an  ordinary  differential  equation  (ODE): 

d£i  =  D  Cl+ 1  -  2C,  +  C,_, 


dt 


Ax2 

,  /^/+ 1  i 

+  I  ^  . - w 

2Ax 


(1  -  <rf)C/+1  +  2 aiCj  -  (1  + 


2Ax 

+  tt/CQ  “  Q  +  2R. 


(3) 


2  Ax 


w 


bb  -  I  (  ~  D*-'  \ 

bb‘  Ax2  +  PA~ ITx - W]~a‘ 

CC.  =  J2-  +  (l  ~  (£i±i  ~  D/-1  _  \ 

Ax2  (  2Ax  A  2Ax  w) 

dd;  =  a,C o  +  ^Ri 


(6a) 

(6b) 

(6c) 

(6d) 


If  we  use  the  known  concentrations  at  WSI  (i.e.,  x  =  0,  first 
node),  C0,  and  at  the  bottom  of  modeled  sediment  column 
(i.e.,  x  =  Lb ,  nth  node),  Ch  as  the  boundary  conditions,  then 

Q- 1  =  c  0;  and  C^„  «  Cb  (7) 

As  a  consequence,  the  series  of  equations  can  be  written  as 
(bb  t  ccx 
aa2  bb2  cc2 


aa^  bbj^i 

aa,  bb, 


CCf 

bbi+i 


cc. 


aan- 1  bb„_  1  cc„_, 

aan  bbn 


X 


'  c, 2 

r  ddl  2 

c2 

dd2 

c„, 

ddj-\ 

c. 

4- 

dd, 

Ci+i 

dd1+  j 

ddn~  i 

k  j 

k  ddn  ; 

=  0  or 


(8a) 


ABC  XC  +  DD  =  0 


(8b) 


where 


in  which  C  and  DD  are  vectors  and  ABC  is  a  tridiagonal 
matrix. 
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Table  2.  Rate  expressions  and  net  reaction  rates  considered  in  this  study 


Notation 


Reaction  rate  expressions: 


R,  =  kG 


[02 


R,  =  kG 


K0l  +  [OJ 
[NO3-] 


K'o,. 


Kno,  +  [NO,“  ]  Kq2  +  [02] 


-  ^NC 


R,  =  kG 


K'. 


[Fe(OH)3]0^r'0i  +  [NO3-]  K'0i  +  [02] 


R<  =  kG- 


[SOj 


K'02 


Kso,  +  [son  J  +  ^Lon,;[Fe(OH)3]  Kf,  +  [N03  ]  K'0l  +  [02] 

[Fe(OH)3]0 

Rs  =  A3[Fe2+][02]  R,  =  *6[NH;][0:]  R7  =  *7[TS][02] 

*8  =  *s[TS][Fe(OH)3]  R9  =  *9[FeS][02]  RUI  =  £10[Fe2+][TS] 

Wand  K°  iS/e  r:^Wgar,lCtrCTC^t,0n  (mo1  C  dm"3-S0,id  s_l)*  K«o»  and  ^so<  are  the  Monod  saturation  constants, 
ana  aD2,  Xno„  and  AFe(0H),  are  the  Monod  inhibition  constants.  ’ 

Equilibrium  expression: 


Net  production  rates: 


[NH4  (ads)]  ^  __  ^^[NH^]  where  p  is  the  sediment  density. 


~  \&s  ~2R6-  207  -  2- - ~Ra 


9 


2*Nor  =  -0.8/?2  +  R 


^Fc(OH),  “  4/?3  +  ^5  “  ^8 


Xi?Fe2,  -  4- 


<P 

1  —  CP  1  —  (p 

- ~R,  + - 

<P  9 


X/?NH4*  —  “(/?,  +  i?2  +  /?3  +  i?4)  —  i?6 


1  ~ 
2  (p 


^R7S  ~  n  7^*4  “  ^7  “  - - “^8  ~  RH 

Z  <p 


^^FeS  ^9  Z  ~  ^10 

I  -  <p 

2*ico,  =  +  7?2  +  i?3  +  7?4) 


7/jvewe  modeling  strategy’  for  rate  determination— The 
objective  of  inverse  modeling  in  the  present  study  was  to 
determine  57?  (and,  thus,  kG)  as  a  function  of  depth  (x).  As 
shown  in  Table  2,  57?  values  are  the  functions  of  kG,  if  the 
values  of  C,  k,  ~  kin,  f,  and  <p  are  known  a  priori.  One 
important  attribute  of  the  study  is  the  a  priori  derivation  of 
Db  and  a:  Whereas  most  previous  studies  treated  DB  and  a 
as  adjustable,  we  assigned  them  a  range  of  a  priori  values 


based  on  the  burrow  and  02  microprofile  data,  as  described 
below.  A  simple  matrix  manipulation  of  (8-2)  yields 

C  =  -ABC-'  X  DD  (9) 

which  indicates  that  the  concentration  of  the  given  species 
spl  at  the  rth  node  can  be  expressed  as  a  function  of  a  priori 
parameters  and  net  reaction  rate,  which  is,  in  turn,  a  function 
of  a  priori  parameters  and  kG, 
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Q5'c ./  =  A(D?-\ „  £>f,  Dg.1,,  af1,  XRf)  (10) 

=  /#g,  cjsl»u.  c£Lre<u,  casL*,.,, . . .)  (ii> 

When  Eqs.  10  and  1 1  are  combined,  kG  is  the  only  un¬ 
known.  It  should  be  noted  that  %Rf'  is  a  function  not  only 
°f  QmLurcd,;  but  also  of  the  measured  concentrations  of  other 
species. 

In  practice,  the  inverse  determination  of  the  kG  profile 
progresses  as  follows.  First,  a  randomly  generated  set  of  ini¬ 
tial  values  are  assigned  to  kGlf  kG2i  .  .  kG„.  It  is  assumed 
that  kG  value  remains  constant  within  each  0.01  m  vertical 
section.  Then,  the  kG  values  are  used  in  Eqs.  9-11  to  cal¬ 
culate  the  concentration  profiles  for  N03 ,  Fe(OH)3,  SOj, 
and  NH4+.  Next,  to  pursue  the  simultaneous  agreement  be- 
tween  Cca!c  and  Cmcilsurcd  for  the  multiple  major  redox  species, 
the  differences  between  Ccale  and  Cmeasured  profiles  are  evalu¬ 
ated  as 


+  y  ISfe*  ~ 

£  cag** 


+  2 


.vLrCJ!! 

"  r'NHf 

where  Qp"ics  denotes  the  highest  measured  concentration  of 
the  given  species  within  the  profile  under  consideration.  Sub¬ 
sequently,  the  kG  values  are  adjusted  to  find  the  optimized 
set  of  values  until  the  value  of  XF  is  minimized  according 
to  the  trust  region  methods  for  nonlinear  minimization  (Co¬ 
leman  and  Li  1996),  using  a  routine  supplied  in  MATLAB 
Optimization  Toolbox  (The  Mathworks,  Inc.).  To  ensure  that 
the  optimization  leads  to  the  global  minima  rather  than  to 
local  minima,  the  process  is  repeated  up  to  50  times,  each 
with  a  new,  randomly  generated  initial  kG  value. 


Forward  calculations  of  depth- concentration  profiles — 
Once  the  depth  profile  of  kG  is  determined,  it  is  then  used 
in  the  forward  modeling  to  calculate  the  depth  profiles  of 
species  listed  in  Table  1 .  In  the  forward  modeling,  Eq.  3  is 
written  for  each  species  using  the  same  Dit  w,  and  a,  values 
used  in  the  inverse  modeling  and  solved  for  C,  by  integrating 
for  time  using  a  stiff  ODE  solver  in  Matlab  (ode23tb). 
is  a  function  of  kG,  which  is  now  a  “known”  parameter, 
and  of  C,  values  of  multiple  species  that  are  now  treated  as 
unknowns.  Because  the  net  rate  of  a  given  species  is  related 
to  the  concentrations  of  multiple  species  (Eq.  1 1),  the  inte¬ 
gration  of  Eq.  3  for  all  species  must  occur  simultaneously 
with  the  rate  expressions  being  coupled  for  each  time  step 
(Boudreau  19966). 


adjacent  to  the  Skidaway  Institute  of  Oceanography,  Geor¬ 
gia,  during  August  2001.  The  sampling  took  place  during 
low  tide  in  the  mid-  to  late-morning  hours.  The  sediment 
temperatures  were  27-28°C,  and  salinity  of  water  from  the 
drainage  creek  was  30.  Data  obtained  from  three  stations  are 
discussed  here,  which  differed  in  the  abundance  of  Ucapug- 
nax  as  well  as  the  density  of  S.  alterniflora  vegetation.  The 
first  station,  designated  BUC,  is  heavily  bioturbated  by  U \ 
pugnax  and  polychaetes  but  contains  no  vegetation.  The  sec¬ 
ond  station,  NUC,  hosts  neither  Spartina  nor  Uca  but  is  pop¬ 
ulated  with  a  small  number  of  polychaetes.  The  third  station, 
NVC,  hosts  abundant  S.  alterniflora  and  a  small  number  of 
polychaetes  but  no  Uca.  The  station  designation  follows  the 
convention  used  by  the  previous  studies  at  SERF  (Kostka  et 
al.  2002a;  Gribsholt  et  al.  2003). 

Ecological  measurements — The  radius,  opening  density, 
and  depth  dependency  of  polychaete  burrows  were  deter¬ 
mined  using  X-radiographs,  as  described  elsewhere  (Furu- 
kawa  et  al.  2001).  The  corer  used  for  X-radiography  was  17 
cm  wide,  2.2  cm  thick,  and  —50  cm  tall.  Because  the  number 
of  Uca  burrows  captured  in  the  corers  was  too  small  for 
statistical  considerations,  depth  dependency  of  Uca  burrows 
was  obtained  from  the  literature  (Basan  and  Frey  1977).  The 
opening  density  and  radii  of  Uca  burrows  were  determined 
directly  in  the  field,  as  described  elsewhere  (Gribsholt  et  al. 
2003)  using  a  25  X  25  cm2  quadrant. 

Porewater  and  solid-phase  geochemistry — The  depth  pro¬ 
files  of  N03-,  S04~,  NH+,  total  salinity  (TS),  Fe2+,  FeS,  and 
Fe(OH)3  were  determined  from  1-  to  2-cm  sections  of  trip¬ 
licate  cores  for  each  station,  as  described  elsewhere  (Kostka 
et  al.  2002a;  Gribsholt  et  al.  2003).  The  02  microprofiles  at 
the  WSI  were  obtained  with  Clark-type  micro  electrodes 
from  freshly  sampled  cores.  The  radial  02  microprofiles  at 
Uca  burrow  walls  were  obtained  by  inserting  the  microelec¬ 
trode  at  a  45°  angle.  The  microelectrodes  were  calibrated  by 
oxygenated  tap  water  whose  02  concentration  was  known 
from  Winkler  titrations  and  by  assuming  that  the  constant 
current  reading  obtained  by  inserting  the  microelectrode  tips 
deeply  into  the  sediments  away  from  burrow  walls  corre¬ 
sponds  to  [02]  =  0  mol  L*1.  It  should  be  noted  that  the 
calibration  was  done  with  the  low  salinity  tap  water;  thus, 
the  data  may  overestimate  the  actual  02  concentrations  in 
pore  water  by  several  percent  (Glud  et  al.  2000). 

Radionuclide  profiles — The  depth  profiles  of  137Cs  were 
obtained  from  centimeter-thick  sections  of  duplicate  cores 
(Alexander  et  al.  1991;  Sommeifreld  et  al.  1999).  Each  sam¬ 
ple  of  dried,  ground  sediment  was  counted  for  -24  h  on  a 
planar  ORTEC  Lo-AX  detector,  to  quantify  ,37Cs  activity  at 
661  KeV.  Self-adsorption  corrections  were  made  using  es¬ 
tablished  techniques  (Cutshall  et  al.  1983).  Errors  were  es¬ 
timated  from  the  1-SD  confidence  level  of  counting  statistics 
to  be  —10%. 


Study  area  and  methods 

We  used  extensive  field  data  obtained  from  muddy  sedi¬ 
ments  of  the  Saltmarsh  Ecosystem  Research  Facility  (SERF) 


Sulfate  reduction  rates  (SSRs) — SRRs  were  determined  on 
10- cm  triplicate  cores  (2  cm  inner  diameter)  with  35S02~ 
(Fossing  and  Jorgensen  1989).  Cores  were  stored  at  the  in 
situ  temperature  for  <2  h  before  injection.  The  35SO^"  was 
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Table  3.  Model  parameters. 


Parameter  Values 


Concentration  parameters 


[NOJ. 

2.5  X  10“6 

[NF;], 

0 

ITS]*  [Fe^J, 

0 

[son„ 

0.029  (BUC); 

[Fe(OH)3]„ 

0.023  (NUC); 
0.023  (NVC) 
0.7 

Transport  parameters 

9 

0.860  (BUC); 

OL  (m2  a*') 

0.821  (NUC); 
0.859  (NVC) 
2.36  X  10-’ 

Dmh  (mJ  s  ') 

2.00  X  109 

(m2  s-) 

2.07  X  10-’ 

D°s o.  (m2  s'1) 

1.11  X  10-’ 

s") 

1.83  X  10-9 

Ohs  (m2  s-‘) 

1.78  X  10-9 

f^Fc(in  (m2  S_I) 

7.36  X  10-'° 

^h2coj  (m2  s  ') 

1.87  X  10-’ 

^FICOj  (nF  S’*) 

1.25  X  10~9 

(m2  s-') 

0.97  X  10“9 

Reaction  parameters 

ks  (L  mol-1  s_l) 

5.9 

k6  (L  mol'1  s“l) 

0.63 

k-j  (L  moP1  s_1) 

0.039 

kx  (L  mol-1  s'1) 

7.0  X  10-4 

k9  (L  mol-1  s~l) 

0.7 

k\a 

0.47 

K0l 

8.0  X  10"6 

•^NCV 

3.0  X  10“5 

Kso, 

0.001 

E Fe(OH), 

2.70  X  10-' 

k. 

1.4 

_ _  Remarks 

Measured  during  1987-1996  at  Skidaway  River  (Verity  2002)  (mol  L~») 
Estimated  (mol  L“l) 

Estimated  (mol  L-1) 

Measured  (this  study)  (mol  L1) 


Estimated  from  Fe(OH)3  and  FeS  profiles  (mol  dm-'3) 
Measured  porosity  (depth  independent) 


Diffusion  coefficients  for  infinite  dilution  is  determined  using  the  tables 
provided  by  Boudreau  (1996a)  with  measured  temperature  (28°C)  and  sa¬ 
linity  (30).  The  tortuosity  correction  follows  the  equations 
Ds  =  D°/P 
P  =  1  -  2  ln(^) 

(Boudreau  199 6a) 


(Zhang  and  Millero  1993) 

(Berg  et  al.  2003) 

(Millero  et  al.  1987) 

(Van  Cappellen  and  Wang  1996;  Boudreau  et  al.  1998) 
(Van  Cappellen  and  Wang  1996;  Boudreau  et  al.  1998) 
(Berg  et  al.  2003) 

Monod  constants  (Boudreau  et  al.  1998) 

Monod  constants  (Boudreau  et  al.  1998) 

Monod  constants  (Boudreau  et  al.  1998) 

Monod  constants  (Boudreau  et  al.  1998) 

NHJ  adsorption  coefficient  (Mackin  and  Aller  1984) 


injected  at  1-cm  intervals,  and  cores  were  incubated  at  in 
situ  temperature  for  3  h.  The  cores  were  then  sectioned  into 
1-2  cm  intervals,  fixed  in  10  ml  20%  ZnAc,  and  frozen. 
Total  reduced  sulfur  (TRS;  i.e.,  acid-volatile  sulfides,  FeS  + 
TS  plus  chromium-reducible  sulfur,  S°  +  FeS2)  was  extracted 
by  a  45 -min  single-step  distillation  with  cold  HC1  (2  mol 
L-1)  and  boiling  Cr2+  (0.5  mol  L~l)  solutions.  Color  obser¬ 
vation  ensured  the  complete  reduction  by  Cr.  TRS  in  the 
distillate  was  determined  spectrophotometrically,  whereas 
TR35S  activity  was  determined  by  a  scintillation  counter. 
Pore  water  was  removed  before  distillations  and  was  also 
analyzed  for  35S02'\  Total  SO2"  was  analyzed  by  ion  chro¬ 
matograph  for  use  in  calculating  rates.  SRRs  were  deter¬ 
mined  as  the  amount  of 35 SO2”  reduced  per  hour,  under  the 
assumption  that  rates  were  constant  during  the  3-h  incuba¬ 
tion. 

Parameter  determination 

Parameters  required  for  the  inverse  and  forward  modeling 
are  listed  in  Table  3.  The  values  in  Table  3  that  are  in  bold¬ 
face  were  determined  through  direct  observations  during  the 


field  study.  The  WSI  is  well  exposed  to  air  or  oxygenated 
estuarine  water;  thus,  the  concentrations  of  reduced  species 
at  the  WSI  were  estimated  to  be  zero.  The  parameterization 
of  Eqs.  9-11  for  N03",  Fe(OH)3,  SO2-,  and  NH+  required 
^measured  profiles  for  02,  N03-,  Fe(OH)3,  SO^,  Fe2+,  TS, 
NHJ,  and  FeS  (see  Table  2).  The  Cmmt  for  02  was  deter¬ 
mined  using  the  directly  measured  microprofiles  plus  burrow 
distribution  statistics.  Cmcasured  for  other  species  was  taken 
from  the  sample  analysis  and  interpolated.  The  <p  was  mea¬ 
sured  in  cores  by  weighing  a  known  volume  of  sediment 
subsamples  before  and  after  drying.  Ds,  which  was  assumed 
to  be  depth  independent,  was  determined  from  its  empirical 
dependence  on  <p  (Boudreau  1996a).  DB  and  a  were  deter¬ 
mined  using  the  burrow  statistics,  and  w  was  determined 
using  radiochemistry. 

The  rate  and  thermodynamic  constants  were  taken  from 
previous  work.  The  values  of  ks  and  k-,  were  calculated  from 
their  empirical  dependence  on  temperature  and  salinity  (Mil¬ 
lero  et  al.  1987;  Zhang  and  Millero  1993).  The  k6  values 
ranging  from  0.16  to  0.63  L  mol-'  s“‘  have  been  used  for 
various  coastal  environments  (Van  Cappellen  and  Wang 
1996;  Berg  et  al.  2003).  The  ks  values  in  previous  work 
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l37Gs  activity  (dpm  g'1) 


Fig.  1.  Depth  profiles  of  excess  l37Cs  activities  at  BUC,  NUC, 
and  NVC.  The  l37Cs  profiles  exhibited  very  little  depth-dependent 
decreases  for  the  sampling  depths  of  the  study.  For  NVC,  the  mea¬ 
sured  depth  profile  is  plotted  together  with  the  calculated  profiles 
according  to  Eq.  (14).  (a)  The  calculated  profile  with  w  =  10  cm 
yr_l.  (b)  The  calculated  profile  with  w  =  3  cm  yr1.  (c)  The  cal¬ 
culated  profile  with  w  =  2  cm  yr A  significant  portion  of  the 
calculated  profile  for  w  =  2  cm  yr’1  falls  outside  of  the  margin  of 
error. 


range  from  1.3  X  lO’4  to  7.0  X  10-"  L  mol"1  S-\  and  k9 
ranges  from  6.0  X  10“3  to  0.7  L  mol'1  s"1  (Van  Cappellen 
and  Wang  1996;  Boudreau  et  al.  1998;  Berg  et  al.  2003).  In 
the  present  study,  the  k6,  ks,  and  k9  values  were  taken  from 
the  highest  values  within  these  ranges,  because  the  high  tem¬ 
perature  of  the  study  site  favored  high  rates.  The  pH-inde- 
pendent  bimolecular  rate  expression  for  FeS  precipitation, 
*io.  has  only  been  tested  recently  in  the  arctic  environment 
(Berg  et  al.  2003).  A  kl0  value  of  0.47  L  mol-1  s_1  was  used 
in  this  study.  A  higher  value  (4.7)  was  also  tested  but  was 
found  to  have  little  effect  on  the  modeled  kG,  Ccalc,  and  TEA 
profiles.  The  value  of  kN  was  taken  from  previous  studies  of 
costal  fine-grained  sediments  (Mackin  and  Aller  1984).  The 
carbon :  nitrogen  ratio  of  labile  OM  was  taken  from  a  recent 
study  conducted  at  SERF  (Gribsholt  et  al.  2003).  The  Mo  nod 
constants  were  taken  from  previous  work  (Boudreau  et  al 
1998). 


Radionuclide  profiles  and  sediment  accumulation  rates — 
Activity  profiles  of  l37Cs  (half-life,  Tos  =  30.2  yr)  from  all 
stations  lacked  depth  dependency  (Fig.  1),  which  could  be 
due  to  rapid  sedimentation  or  intense  bioturbation.  Mean¬ 
while,  X-radiographs  exhibited  visible  horizontal  units, 
which  indicates  that  bioturbation  has  not  been  intense 
enough  to  completely  destroy  sedimentary  structures.  Thus, 
the  lack  of  depth  dependency  must  be  due  to  the  rapid  sed¬ 
imentation.  By  applying  the  steady-state  diffusion-decay 
equation, 


w— - AC  =  0 

dx 


(13) 


where  A  is  the  decay  constant  (0.693/7’O5)  to  the  l37Cs  profile 
from  NVC,  under  the  assumption  that  D„  =  0,  one  obtains 
the  l37Cs  concentration  as  the  function  of  depth, 


[l37Cs]  =  [137Cs],=0exp^-^  (14) 

By  applying  the  w  values  of  2,  3,  and  10  cm  yr-1  in  Eq.  14, 
the  modeled  profiles  (Fig.  1;  NVC)  are  obtained.  The 
[137Cs]x  =  0  value  is  adjusted  for  each  case  to  seek  a  fit  to 
the  observed  profile.  The  w  value  of  10  cm  yr-1  yields  a 
profile  that  is  nearly  independent  of  depth  and  fits  most  of 
the  observed  data  points.  Whereas  w  =  3  cm  yr-1  yields  a 
profile  with  some  depth  dependency,  the  modeled  profile  still 
falls  within  the  counting  errors  of  most  of  the  data  points. 
However,  w  =  2  cm  yr-1  deviates  further  from  the  data 
points  with  more  significant  depth  dependency.  As  a  con¬ 
sequence,  we  used  w  values  of  10  and  3  cm  yr"1.  The  effect 
of  higher  (i.e.,  20  cm  yr-1)  w  value  is  also  examined. 


Burrow  parameters  and  irrigation  coefficient,  a — The  ir¬ 
rigation  coefficient  was  estimated  from  the  burrow  obser¬ 
vations  using  existing  formulae  (Koretsky  et  al.  2002)  in 
which  ecological  measurements  such  as  the  number  and  size 
of  burrows  can  be  related  to  a  as 


(15) 


Sv  is  the  total  burrow  surface  area  divided  by  the  total  vol¬ 
ume  of  sediments  modeled  for  the  given  depth  interval,  and 
L  is  the  radial  diffusion  length  scale  which  is  the  difference 
between  r,  the  radial  distance  from  the  burrow  axis  at  which 
the  solute  concentration  reaches  the  radially  averaged  con¬ 
centration  for  the  bulk  sediment  at  a  given  depth,  and  r„  the 
equivalent  burrow  radius  (Aller  1980;  Boudreau  1984).  The 
equivalent  burrow  radius  is  related  to  the  total  burrow  sur¬ 
face  area  and  density  at  a  given  depth  in  the  equation 


rt(x) 


AS, 

Ax 

2  rr«to,(x) 


(16) 


where  nM(x)  is  the  total  number  of  burrows  found  at  the 
given  depth  x  within  the  unit  area. 

The  values  of  L,  S„  and  r„  as  well  as  their  depth 
dependency,  are  shown  in  Table  4  for  all  stations.  The  poly- 
chaete’s  contributions  to  5j,  and  «10t  were  directly  determined 
from  the  X-radiography  data  of  burrow  numbers,  size,  and 
tilt  angles  (Furukawa  et  al.  2001),  whereas  the  contributions 
from  U.  pugnax  were  estimated  by  assuming  that  all  Uca 
burrows  observed  at  the  WSI  inside  the  25  X  25  cm2  quad¬ 
rant  had  1  depth  extent  of  0.2  m  and  a  radius  of  0.0075  m, 
which  are  the  average  depth  (0.2  ±  0.025  m)  and  radius 
(0.0075  ±  0.0025  m)  observed  by  previous  studies  (Teal 
1958;  Allen  and  Curran  1974;  Basan  and  Frey  1977).  The 
total  Sv  was  simply  the  sum  of  polychaete  and  Uca  contri¬ 
butions. 

The  value  of  L  for  Uca  burrows  was  obtained  by  (1)  as¬ 
suming  the  radial  02  penetration  distance  to  be  equal  to  L 
(Koretsky  et  al.  2002)  and  (2)  averaging  the  direct  measure¬ 
ments  of  radial  O,  penetration  depths  obtained  along  Uca 
burrows  at  the  study  site  (0.82  X  lO"3  [±3.3  X  10'4]  m). 
The  radial  02  penetration  depth  along  polychaete  burrows 
was  not  directly  measured.  The  L  (F-  r,)  value  was  instead 
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determined  from  the  following  ODE,  boundary  conditions, 
and  an  additional  equation  that  assumed  that  Ds  and  the  02 
consumption  rate  (R)  are  identical  between  the  sediments 
surrounding  Uca  burrows  and  surrounding  polychaete  bur¬ 
rows: 


l_d 
r  dr 


+  R  =  0, 


0, 


dC 

dr 


=  0  (17) 


(Cai  and  Sayles  1 996).  For  BUC,  a  was  calculated  separately 
for  Uca  and  polychaetes,  according  to  Eq.  15,  before  being 
added  together  to  obtain  the  total  a. 

It  is  expected  that  phytoirrigation  is  important  in  the  pres¬ 
ence  of  S.  alterniflora .  However,  the  intensity  of  irrigation 
cannot  be  quantified  using  available  data,  because  no  infor¬ 
mation  is  currently  available  that  links  the  intensity  of  phy¬ 
toirrigation  to  the  ecological  data  related  to  vegetation.  As 
a  consequence,  a  for  NVC  was  derived  using  the  burrow 
parameters  in  the  same  manner  as  for  BUC  and  NUC,  under 
the  assumption  that  phytoirrigation  is  insignificant  compared 
with  the  macrofaunal  irrigation.  The  consequence  of  this  as¬ 
sumption  is  discussed  later. 


Biodiffusion  coefficient ,  DB— We  obtained  estimated  DB 
values  by  assuming  that  the  majority  of  particle  displace¬ 
ment  activity  is  due  to  the  particle  ingestion/egestion  by 
macrofauna.  The  magnitude  of  particle  ingestion/egestion 
has  been  found  to  be  roughly  proportional  to  the  size  of 
animals  (Cammen  1980).  The  volumetric  size  of  animals  can 
be  regarded  as  proportional  to  the  cube  of  the  burrow  radius. 
In  turn,  the  unit  area-integrated  burrow  radius  can  be  re¬ 
garded  as  proportional  to  the  cube  of  burrow  surface  area 
summed  over  all  burrows  divided  by  the  total  volume  of 
sediments  modeled  within  the  given  depth  interval: 

db  -  y&l 

The  value  of  scaling  factor,  y,  is  adjusted  so  that  DB  averaged 
for  the  top  10  cm  of  heavily  bioturbated  BUC  falls  within 
the  empirically  compiled  DB  range  of  10  <  DB  <  1,000  cm2 
yr~!  for  bioturbated  and  rapidly  accumulating  sediments 
(Boudreau  1994).  In  practice,  this  meant  that  the  value  of  y 
was  varied  between  1.8  X  10~14  and  1.8  X  10“52. 


Depth-0 2  concentration  profiles— The  Cmeasurcd  profiles  for 
02  were  obtained  by  combining  measured  vertical  micropro¬ 
files  at  WSI  and  calculated  effects  of  radial  02  penetrations 
along  burrow  walls.  The  effect  of  radial  penetration  on  the 
vertical  profiles  was  determined  by  radially  integrating  the 
penetration  profile  perpendicular  to  the  Uca  burrows  and  av¬ 
eraging  over  the  unit  volume  of  sediments: 


C  -  N 


2  7 rC(r)  dr 


(19) 


N  is  the  number  of_burrows  per  unit  lateral  cross  section 
area  of  sediments,  r]  is  the  average  burrow  radius  for  the 
given  depth,  and  r,  +  L  is  the  average  02  penetration  radius. 


C(r)  is  the  function  for  the  average  radial  02  concentration 
profile  perpendicular  to  a  burrow  wall  away  from  the  WSI 
and  is  derived  by  solving  Eq.  17  for  C.  The  net  vertical  02 
profiles  were  determined  by  combining  the  02  from  vertical 
penetration,  radial  02  penetration  along  Uca  burrows,  and 
radial  02  penetration  along  polychaete  burrows.  The  poly¬ 
chaete  burrow  contribution  was  reduced  to  20%  of  the  values 
derived  from  Eq.  19,  because  polychaete  irrigation  is  inter¬ 
mittent  (Kristensen  2001).  Uca  burrows,  with  their  large 
openings  and  simple  geometry,  were  assumed  to  be  fully 
oxygenated  by  bioirrigation  and  physical  pumping. 

Results 

Inverse  modeling  yielded  the  kG  profiles  shown  in  Fig.  2. 
For  each  station,  different  DB  (i.e.,  different  y  in  Eq.  18;  y 
=  1-8  X  10"14,  1.8  X  10~13,  and  1.8  X  10~12)  and  w  (i.e.,  w 
“  3,  10,  and  20  cm  yr"1)  resulted  in  kG  profiles  that  are 
similar.  On  the  other  hand,  the  comparison  between  stations 
reveals  considerable  differences.  The  OM  degradation  rate 
is  modeled  to  be  significantly  higher  at  BUC  and  NVC  than 
at  NUC.  The  results  also  show  that  a  major  portion  of  OM 
at  NUC  is  remineralized  within  —1  cm  of  the  WSI,  whereas 
the  remineralization  occurs  well  below  WSI  at  BUC  and 
NVC,  regardless  of  the  parameter  values  chosen. 

The  inversely  determined  kG  profiles  were  subsequently 
used  in  the  forward  modeling  to  calculate  the  depth-concen¬ 
tration  profiles  of  major  redox  species  and  to  confirm  the 
agreement  between  measured  and  modeled  profiles  (Fig.  3). 
At  BUC,  the  kG  profiles  calculated  using  the  wide  range  of 
w (3,  10,  and  20  cm  yr_l)  yielded  very  similar  Ccalc  profiles, 
whereas  the  choice  of  y  affected  the  results.  The  results  show 
that  y  =  1.8  X  10-12,  which  yields  an  average  DB  for  the 
top  10  cm  at  BUC  of  1,000  cm2  yr-1,  generates  good  agree¬ 
ment,  whereas  the  smaller  y  (i.e.,  smaller  DB)  resulted  in  the 
modeled  Fe-(oxy)hydroxide  concentrations  to  reach  zero 
within  upper  ~ 5  cm  of  sediments.  The  measurements  indi¬ 
cate  a  significant  presence  of  Fe-(oxy)hydroxides  throughout 
the  sampling  depths.  In  NUC,  the  choice  of  w  affected  the 
calculated  depth  profiles  whereas  the  choice  of  y  had  neg¬ 
ligible  effect.  Model  calculations  yielded  reasonable  agree¬ 
ment  with  the  measured  profiles  when  the  low  sedimentation 
rate  (w  =  3  cm  yr~J)  was  used  in  the  calculations.  Only  the 
results  using  the  kG  profiles  calculated  with  w  =  3  cm  yr 1 
and  y  =  1.8  X  10“12  are  shown  for  NVC. 

The  forward  model  was  also  used  to  determine  the  depth 
profiles  of  OM  degradation  rate  by  each  TEA,  Rlt  R2. ,  R3, 
and  R4 ,  according  to  the  formulae  in  Table  2.  The  inversely 
determined  kG  profiles  were  used  with  the  forward  model- 
determined  TEA  concentrations  in  the  formulae.  Once  these 
were  determined,  the  rates  of  SOf  and  Fe(III)- 
(oxy)hydroxide  consumption  due  to  OM  remineralization 
could  be  calculated  as  follows  according  to  the  reaction  stoi¬ 
chiometry 

1  —  (p  1 

=  — “ — ^  (mo'  L-porewater-1  s~')  (20) 

■KfcTohij  =  (mol  drrr3-so!ids  s'1)  (21) 

The  and  profiles  are  plotted  in  Fig.  4,  together 
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Table  4.  Burrow  parameter  summary. 

BUC 


Polychaete  Total 


Depth 

(m) 

Density 

WU» 

(m“2) 

Radius 

^lUca 

(m) 

Burrow 

surface 

area 

DSv/Dx 

(m) 

Diffusion 

length 

scale 

L 

(m) 

Density 

*V>ly 

(nr2) 

Radius 

ripoly 

(m) 

Burrow 

surface 

area 

DSvIDx 

(m) 

Diffusion 

length 

scale 

L 

(m) 

Burrow 

surface 

area 

DSv/Dx 

(m) 

(DSv/L)Un 

+ 

(DSv/LU 

0.00 

131 

0.0075 

6.17 

7.00X10"4 

713 

5.58X10"4 

2.50 

5.05X10"4 

8.67 

13,763 

0.01 

131 

0.0075 

6.17 

7.00X10-" 

713 

5.58X10-4 

2.50 

5.05X10-4 

8.67 

13,763 

0.02 

131 

0.0075 

6.17 

7.00X10"4 

446 

7.89X10-4 

2.21 

5.45  XI  O'4 

8.38 

12,866 

0.03 

131 

0.0075 

6.17 

7.00X10-4 

357 

5.62X10-4 

1.26 

5.06xi0"4 

7.43 

11^303 

0.04 

131 

0.0075 

6.17 

7.00X10”4 

178 

9.04X10-4 

1.01 

5.60X  10~4 

7.18 

10*615 

0.05 

131 

0.0075 

6.17 

7.00X10-4 

267 

1.01X10-3 

1.69 

5.71  X  10~4 

7.86 

1 1 ,769 

0.06 

131 

0.0075 

6.17 

7.00X10-" 

891 

1.05X10"3 

5.85 

5.75  XIO'4 

12.02 

18,986 

0.07 

131 

0.0075 

6.17 

7.00X10-4 

980 

8.56X10-4 

5.27 

5.54X10-4 

11.44 

18*325 

0.08 

131 

0.0075 

6.17 

7.00X10-4 

446 

8.03X10'4 

2.25 

5.47X10”4 

8.42 

12,925 

0.09 

131 

0.0075 

6.17 

7.00X10-4 

535 

1.07  XIO-3 

3.59 

5.77X1 0"4 

9.76 

15,031 

0.10 

131 

0.0075 

6.17 

7.00X10  4 

178 

7.60X1 0-4 

0.85 

5.41  X10  4 

7.02 

10,382 

0.11 

131 

0.0075 

6.17 

7.00X10-4 

89 

5.55X10-4 

0.31 

5.04X1 0“4 

6.48 

9,426 

0.12 

131 

0.0075 

6.17 

7.00X10"4 

624 

8.70X10"4 

3.41 

5.56X10"4 

9.58 

14,948 

0.13 

131 

0.0075 

6.17 

7.00X10-4 

267 

6.80X10'4 

1.14 

5.28  X10"4 

7.31 

10,970 

0.14 

131 

0.0075 

6.17 

7.00X10-4 

267 

5.96X10-4 

1.00 

5.13X10'4 

7.17 

10,761 

0.15 

131 

0.0075 

6.17 

7.00X10-4 

357 

9.77X10”4 

2.19 

5.68X10-4 

8.36 

12,667 

0.16 

131 

0.0075 

6.17 

7.00X1 0-4 

535 

1.89X10-3 

6.35 

6.28  X10“4 

12.52 

18,914 

0.17 

131 

0.0075 

6.17 

7.00X104 

178 

7.96X10  4 

0.89 

5.46X10  4 

7.06 

10,441 

0.18 

131 

0.0075 

6.17 

7.00X1 0-4 

178 

6.71X1  O'"4 

0.75 

5.26X1 0-4 

6.92 

10,235 

0.19 

131 

0.0075 

6.17 

7.00X1  O'4 

178 

8.05  X10"4 

0.90 

5.47X10-4 

7.07 

10,455 

0.20 

131 

0.0075 

6.17 

7.00X1 0~4 

178 

5.10X10”4 

0.57 

4.94  XIO-4 

6.74 

9,965 

with  the  experimentally  measured  SRR  values.  The  mea¬ 
sured  SRR  is  expressed  in  terms  of  mol  L  pore  water-1  s“*, 
to  allow  a  direct  comparison  with  the  model-calculated 
^scm-  For  BUC,  the  measured  SRR  and  model-calculated 
^so  4  are  in  agreement.  Model-calculated  R™1  profiles  would 
be  three  to  four  times  greater  if  the  smaller  y  values  were 
used,  generating  poor  agreement  with  the  measured  SRR. 
The  model-calculated  profile  is  smaller  than  measured 
SRR  for  NUC.  For  NVC,  the  net  R§0  4  and  SRR  for  the  upper 
9  cm  of  the  sediments  are  in  general  agreement,  although 
the  depth-dependent  fluctuation  found  in  the  model  results 
was  not  observed  in  the  field. 

The  relative  importance  of  Fe(III)-(oxy)hydroxides  as  a 
TEA  was  calculated  to  be  greater  at  the  heavily  bioturbated 
BUC  than  at  NUC  or  NVC  (Fig.  5).  At  NUC  and  NVC, 
model  results  show  that  SO^”  is  responsible  for  the  majority 
of  OM  degradation.  In  contrast,  a  significant  fraction  of  OM 
degradation  proceeds  via  reduction  of  Fe(III)- 
(oxy)hydroxides  at  BUC. 

Discussion 

Effect  of  macrobenthos  on  kG — The  model  results  indi¬ 
cate  that  macrobenthos  affect  the  chemical  mass  transfer  in 
Skidaway  SERF  sediments.  A  comparison  of  BUC  and  NUC 
profiles  (Fig.  2)  illustrates  the  effect  of  bioturbation.  The  net 
rate  of  OM  degradation  is  much  greater  (—8  times)  in  the 
bioturbated  sediments  of  BUC  than  in  NUC.  In  addition, 
most  of  the  OM  degradation  occurs  within  the  immediate 


vicinity  of  WSI  in  NUC,  whereas  it  takes  place  well  into  the 
deeper  part  of  the  sediments  in  BUC.  Although  the  observed 
lack  of  SO}~  and  Fe(OH)3  depletion  at  BUC  (Fig.  3)  may 
suggest  a  lack  of  significant  OM  degradation,  the  modeling, 
with  its  forcing  of  the  high  rate  of  biological  transport,  re¬ 
veals  the  reason  for  the  lack  of  TEA  depletion.  TEAs  pen¬ 
etrate  deeper  at  BUC  not  because  they  are  underutilized  but 
because  they  are  transported  and  regenerated. 

The  comparison  of  kG  profiles  from  NUC  and  NVC  (Fig. 
2)  makes  clear  that  vegetation  also  enhances  the  sedimentary 
OM  degradation  rates:  the  net  OM  degradation  is  greater 
(—14  times)  in  the  heavily  vegetated  sediments  of  NVC  than 
in  NUC.  In  addition,  OM  degradation  takes  place  well  into 
the  deeper  part  of  the  sediments  in  NVC,  whereas  it  is  re¬ 
stricted  to  the  immediate  vicinity  of  WSI  in  NUC. 

Parameter  selections — We  examined  a  range  of  sediment 
accumulation  rates  (w)  and  scaling  factors  for  the  biodiffu¬ 
sion  coefficient  (y).  Appropriate  values  for  these  parameters 
may  be  selected  by  comparing  the  model  results  with  the 
measurements  from  this  and  previous  studies.  A  recent  study 
within  SERF  at  a  bioturbated  station  similar  to  BUC  has 
shown  that,  at  the  depth  of  5  cm,  >50%  of  OM  remineral¬ 
ization  occurs  through  the  reduction  of  Fe(III)  minerals 
(Gribsholt  et  al.  2003).  The  model-calculated  relative  im¬ 
portance  of  TEAs  at  BUC  agrees  with  this  result,  showing 
that  —50%  of  the  OM  degradation  proceeds  through 
Fe(OH)3  reduction  at  5-6  cm  below  WSI  (Fig.  5).  The  lower 
y  values,  y  =  1.8  X  10~14  and  1.8  X  10“13,  would  decrease 
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Table  4.  Extended. 


NUC 


NVC 


Polychaete 


Burrow 

Diffision 

surface 

length 

area 

scale 

Density 

Radius 

DSvIDx 

L 

^poly 

^ipoly 

(m) 

(m) 

201 

5.78  XI 0-4 

0.73 

5.09  XIO"4 

201 

5.78X10-" 

0.73 

5.09X10-4 

334 

9.12X10-4 

1.92 

5.61  X10~4 

67 

6.85  XIO'4 

0.29 

5.29X10-" 

67 

5.00X1 0-* 

0.21 

4.91  X10-" 

134 

5.23X1 0“4 

0.44 

4.97  X10-" 

401 

8.16X10*4 

2.05 

5.49X10-" 

334 

6.35X1 0-" 

1.33 

5.20  XIO"4 

201 

1.57X10-’ 

1.98 

6.13X10-" 

267 

8.28  XIO”4 

1.39 

5.50X10-4 

67 

5.18X10-4 

0.22 

4.96X104 

0 

0.00 

0 

0.00 

0 

0.00 

0 

0.00 

0 

0.00 

0 

0.00 

0 

0.00 

0 

0.00 

134 

7.06X10"4 

0.59 

5.32E-04 

134 

5.03  X10"4 

0.42 

4.92E-04 

Denstiy 

”poiy 

446 

446 

1070 

1693 

1337 

802 

1426 

624 

267 

267 

624 

446 

178 

267 

267 

446 

267 

446 

0 

0 

89 


Polychaete 


Radius 

ripo!y 

Burrow 

surface 

area 

DSvIDx 

(m) 

5.63X10-" 

1.57 

5.63  XIO"4 

1.57 

1.02X10-’ 

6.87 

7.78  X10-4 

8.27 

7.80X10-" 

6.55 

8.24X10-" 

4.15 

1.01  xio-3 

9.09 

1.15X10"3 

4.51 

7.06X10-" 

1.19 

8.26X10-" 

1.39 

9.50X10  4 

3.72 

3.15X10-3 

8.80 

6.07X10-" 

0.68 

7.81  X10-" 

1.31 

7.42X1 0"4 

1.25 

1.08  XIO’3 

3.03 

5.32X10-" 

0.89 

9.82X10  4 

2.75 

0.00 

0.00 

1.39X10-’ 

0.78 

Diffusion 

length 

scale 

L 

(m) 

5.06  X10“4 
5.06X10-4 
5.73X10”4 
5.43  XIO"4 
5.44X1 0~4 
5.50X10’4 
5.72X1 0"4 
5.85X10"4 
5.32X10-4 
5.50X10~4 
5.65X10  4 
6.63  XIO"4 
5.15X10"4 
5.44X1 0”4 
5.38X10'4 
5.78X10-4 
4.99X1 0"4 
5.69X10  4 


6.02X10-4 


the  relative  importance  of  Fe(OH)3  in  this  depth  interval  to 
~0%  and  3%,  respectively.  The  comparison  of  model-cal¬ 
culated  and  measured  SRRs  (Fig.  4)  reveals  that  measured 
SRR  at  BUC  are  reasonably  matched  by  using  y  —  1.8  X 
10-'2.  Moreover,  the  comparison  of  concentration  profiles  at 
BUC  (Fig.  3)  shows  that  the  modeling  with  a  smaller  D„ 
(i.e.,  y  =  1.8  X  10"14  ~  ~13)  would  result  in  the  exhaustion 
of  Fe(III)-(oxy)hydroxides  within  the  upper  ~5  cm  of  sed¬ 
iments,  even  though  significant  amount  is  observed  through¬ 
out  the  sampling  depth.  Thus,  y=  1.8  X  ICM2  is  likely  the 
most  appropriate  choice. 

The  forward-model  simulations  of  concentration  profiles 
using  the  high  sediment  accumulation  rate  values  (i.e.,  10 
and  20  cm  yr-1;  Fig.  3)  result  in  the  deeper  penetration  of 
Fe(III)-(oxy)hydroxides  in  comparison  to  the  directly  mea¬ 
sured  profiles  at  NUC.  Thus,  the  w  values  of  3  cm  yr’1  are 
likely  the  more  appropriate  choice. 

The  combined  use  of  w  =  3  cm  yr-1  and  y  =  1.8  X  10-12 
yields  agreement  between  measured  and  modeled  concentra¬ 
tion  profiles  at  all  stations  (Fig.  3).  There  is  reasonable 
agreement  between  measured  and  model-calculated  SRR  for 
BUC.  The  depth-integrated  measured  SRR  and  modeled 
■^soV  f°r  the  upper  9  cm  of  sediments  are  calculated  as 
follows: 

Jr0.09 

dx  (mol  m~2  s~‘)  (22) 

0 

The  results,  shown  in  Table  5,  are  comparable  for  BUC  and 
NVC. 


Effect  of  C:N ratio  in  NUC— The  NUC  modeling  results 
exhibit  the  underestimation  of  R°qI~  compared  with  the 
measured  SRR  (Fig.  4).  Also,  the  forward  modeling  reveals 
that  the  simultaneous  agreement  between  measured  and 
model-calculated  concentration  profiles  is  not  achieved  for 
NUC  as  S02~  is  overestimated  by  the  model  (Fig.  3).  This 
discrepancy  may  be  due  to  the  model  use  of  an  inappropri¬ 
ately  low  C  :  N  ratio.  A  low  C  :  N  ratio  in  the  NH+  profile 
matching  processes  of  inverse  modeling  would  result  in  the 
underestimation  of  kG,  which  would  lead  to  the  underesti¬ 
mation  of  R™\~.  Given  the  lack  of  fresh  OM  input  from 
vegetation  and  bioturbation,  the  penetration  of  labile,  N-rich 
OM  is  limited  at  NUC  except,  for  the  immediate  vicinity  of 
WSL  The  decomposing  Spartina  plants,  whose  C  :  N  ratio 
has  been  reported  to  be  32-54  at  nearby  Sapelo  Island 
(Haines  et  al.  1977),  would  result  in  a  high  C  :  N  ratio  in 
the  absence  of  more  labile,  N-rich  OM.  By  using  the  high 
end  value  of  C  :  N  =  54  instead  of  12.5  at  NUC  in  con¬ 
junction  with  w  =  3  cm  yr~'  and  y  =  1.8  X  1012,  we 
obtained  forward-model  results  (Fig.  6)  in  which  the  agree¬ 
ment  between  measured  and  model-calculated  profiles  was 
improved.  The  simulation  with  C  :  N  —  54  also  yielded  an 
improved  match  between  measured  SRR  and  model-calcu¬ 
lated  R™ l~  (Fig.  7,  Table  5).  The  depth-integrated  SRR  at 
NUC  using  the  adjusted  C  :  N  were  also  comparable  to  the 
values  measured  in  August  2000  at  a  Ucci-  and  Spartina-free 
part  of  the  Skidaway  SERF  (Kostka  et  al.  2002a). 

The  improved  agreement  suggests  that  the  C  :  N  ratio  of 
degrading  OM  at  NVC  is  likely  to  be  closer  to  54  than  to 
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kG  (X  10-6  mol  dm-3  s*1) 
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Fig.  2.  The  model-calculated  depth  profiles  of  OM  degradation 
rate  (kG)  for  BUC,  NUC,  and  NVC.  In  the  legends,  wl,  w2y  and 
w3  denote  model  runs  using  w  =  3,  10,  and  20  (cm  yr~l),  respec¬ 
tively,  and  y\,  y2}  and  y3  represent  runs  with  y  =  1.8  X  10“14,  1.8 
X  1 0~13,  and  1.8  X  10“12,  respectively.  At  NUC,  the  OM  degra¬ 
dation  is  fastest  in  the  immediate  vicinity  of  WSI  (0-0.01  m)  than 
in  the  deeper  part  of  the  sediments.  On  the  other  hand,  at  BUG  and 
NVC,  the  OM  degradation  occurs  in  the  deeper  parts  of  the  sedi¬ 
ments,  and  the  net  OM  degradation  is  greater  than  NUC. 


12.5.  The  C  :  V  ratio  of  remineralizing  OM  is  different  from 
the  bulk  C  :  N  ratio  in  any  given  sample  (Berg  et  al.  2003), 
and  there  is  no  direct  method  that  can  quantify  the  C  :  N 
ratio  of  degrading  OM  in  situ.  The  parameter  fitting  ap¬ 
proach  taken  here  can  be  a  reasonable  method  when  one 
needs  to  determine  the  C  :  N  ratio  of  degrading  OM. 

Integrated  TEA  utilization  rates  versus  TEA  fluxes— TEA 
utilization  rate  profiles  can  be  used  to  calculate  the  net,  in¬ 
tegrated  TEA  utilization  rates  for  the  upper  9  cm  of  sedi¬ 
ments  at  BUC,  NUC,  and  NVC.  For  Fe(OH)3, 

J'"  0.09 

^Feam  dx  (mol  m-2  s"1)  (23) 

o 

The  value  for  each  station  is  shown  in  Table  5.  The  inte¬ 
grated  value  at  BUC  is  comparable  to  4.0  X  10~6  mol  m~2 
s_l,  the  value  found  previously  by  an  incubation  for  a  bio- 
turbated  and  vegetated  substation  in  SERF  (Kostka  et  al. 
2002a).  The  flux  of  Fe(III)-(oxy)hydroxides  to  the  SERF 
sediments,  FPc(im  mol  m-2  s"1,  can  be  estimated  from  the  total 


Fig.  3.  The  comparison  of  model-calculated  and  measured  con¬ 
centration-depth  profiles  of  major  redox  species.  For  BUC,  different 
y  values  yielded  different  results,  as  shown.  For  NUC,  the  selection 
of  sediment  accumulation  rates  (w)  affected  the  calculation  results, 
as  shown. 


solid-phase  iron,  XFe  (moles  dm~3-solids),  and  sediment  ac¬ 
cumulation  rates,  w  (m  s~J),  as  follows: 

^Fc(iii)  =  103w(l  -  <p)  XFe  (24) 

where  XFe  is  the  total  solid-phase  iron  (=  Fe(III)  +  Fe(II)) 
averaged  for  the  entire  depth  of  sediment  columns  studied. 
The  directly  measured  XFe  are  similar  at  all  three  stations 
(i.e.,  XFe  ~  0.7  moles  L-solids-1).  Consequently,  the  flux  of 
Fe(III)-(oxy)hydroxides  to  either  station  is  the  function  of  w 
and  <p,  and  the  value,  nearly  identical  for  all  stations,  is 
shown  in  Table  5. 

The  comparison  of  XRFe(III)  and  FPem  reveals  that,  at  NUC, 
the  flux  of  Fe(III)-(oxy)hydroxide  due  to  sediment  accu¬ 
mulation  is  closely  matched  by  the  calculated  cumulative 
rate  of  its  consumption  because  of  OM  degradation  within 
upper  9  cm  of  the  sediments.  The  lack  of  Fe  recycling  in 
turn  increases  the  relative  importance  of  SO2”  as  the  TEA. 
A  recent  study  at  a  Uca-  and  Spartina- free  part  of  the  Ski- 
daw ay  SERF  found  that  >70%  of  OM  remineralization  was 
through  the  reduction  of  S02“  (Kostka  et  al.  2002a).  At 
NVC,  the  integrated  Fe(III)  consumption  rate  is  also  similar 
to  the  flux. 

In  contrast,  at  BUC,  the  calculated  cumulative  consump- 
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Fig.  4.  Model-calculated  depth  profiles  of  Fe(III)-(oxy)hydroxide 
reduction  rates  and  SOJ~  reduction  rates  caused  by  OM  reminer¬ 
alization.  The  measured  SRR  profiles  are  also  plotted  for  compar¬ 
ison.  It  should  be  noted  that  SRR  is  expressed  in  terms  of  the  rate 
of  SO Y  consumption  in  the  unit  volume  of  pore  water  ( jamol  L_l 
s~l).  For  BUG,  the  choice  of  smaller  y  values  (y  =  1.8  X  10"14. 
and  1.8  X  10~13)  resulted  in  the  SO^~  reduction  rate  that  is  three 
to  four  times  greater  than  the  measured  SRR. 

tion  rate  of  Fe(III)-(oxy)hydroxides  in  the  upper  9  cm  of 
sediments  is  30  times  greater  than  the  flux,  which  suggests 
that  each  Fe  atom  experiences  cycling  between  reducing  and 
oxidizing  environments  at  least  -30  times  via  bioturbation 
before  it  is  permanently  buried  below  the  zone  of  biological 
mixing.  This  value  is  underestimating  the  magnitude  of  re- 


Relative  TEA  usage 


Fig.  5.  The  relative  importance  of  02,  N03~,  Fe(III),  and  SOJ"  as 
TEAs  in  the  OM  remineralization,  plotted  as  a  function  of  depth. 
These  were  calculated  by  normalizing  the  model-calculated  Ru 
R3,  or  R,  by  (*,  +  R2  +  R2  +  RA). 


cycling,  because  the  zone  of  biological  mixing  at  BUC  is 
>9  cm  (Table  4).  This  significant  recycling  results  in  the 
enhanced  relative  importance  of  Fe(III)  as  a  TEA  at  BUC. 
The  rate  expression  for  Fe(OH)3  production  (Table  2)  indi¬ 
cates  that  it  depends  on  one  reaction,  the  aerobic  oxidation 
of  Fe2+  (i.e.,  R5).  The  high  bioturbation  and  irrigation  rates 
are  the  key  to  maintaining  the  deep  penetration  of  02  that 
allows  the  continuous  production  of  Fe(OH)3.  The  Fe  recy¬ 
cling  of  —30  times  is  approaching  the  value  of  100-300 
determined  in  the  shallow  marine  environment  of  Skagerrak 
off  the  coast  of  Denmark  (Canfield  et  ai.  1993).  Meanwhile, 
in  sediments  with  limited  particle  mixing  (i.e.,  NUC  and 
NVC),  solid-phase  TEAs  such  as  Fe(III)-(oxy)hydroxides 
are  reduced  once  and  are  not  reoxidized  to  participate  in 
repeated  reduction. 

The  model  provides  independent  confirmation  that  biotur¬ 
bation  drives  particle  and  solute  mixing  at  BUC.  Previous 
experimental  studies  (Thamdrup  2000;  Kostka  et  al. 
2002 a,b\  Gribsholt  et  at.  2003)  have  indicated  that  Fe(III) 
respiration  is  sustained  by  metal  oxides  supplied  via  biotur¬ 
bation.  Because  Fe(III)  minerals  becomes  less  reactive  with 
age  after  precipitation,  it  may  also  be  speculated  that  rapid 
recycling  due  to  bioturbation  keeps  Fe(III)  minerals  in  a 
poorly  crystalline  state,  exerting  a  positive  feedback  on  mi¬ 
crobial  Fe(III)  respiration  (Thamdrup  2000).  This  positive 
feedback  would  also  support  the  enhanced  OM  degradation 
rates  under  redox-oscillating  environments  (Aller  1994). 

Effects  of  phytoirrigation— To  examine  whether  phytoir¬ 
rigation  is  quantitatively  significant  for  the  transport  of  major 
solute  species,  additional  model  runs  were  conducted  for 
NVC.  It  should  be  noted  that  the  elevated  a  would  affect  all 
species  equally  and  would  not  simulate  the  uptake  of  specific 
solute  species.  One  simulation  was  conducted  with  the 
depth-independent  a  whose  value  was  derived  by  first  taking 
the  average,  macrofauna-only  a  values  for  each  species  in 
the  upper  10  cm  of  NVC  sediments,  and  increasing  it  by 
three  times.  Another  simulation  was  conducted  by  using  ten 
times  the  average,  macrofauna-only  a  value.  The  model-cal¬ 
culated  profiles  using  the  elevated  a  values  signifi¬ 

cantly  overestimate  the  measured  SRR  values  at  NVC, 
which  indicates  that  the  elevated  a  values  do  not  appropri¬ 
ately  depict  the  irrigation  in  NVC.  This  implies  that  the  net 
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Table  5.  Comparison  of  depth-integrated  measured  SRR  and  model-calculated  as  well  as  the  depth-integrated  Fe(III)  flux  and 
model-calculated  ^fck>h)}  f°r  the  upper  0.09  m  of  sediments.  The  consumption  rates  were  determined  by  the  model  using  the  parameter 
values  of  w  =  3  cm  yr1  and  y  =  1.8  X  10  l2,  whereas  the  flux  was  estimated  from  the  measured  total  iron  concentrations  and  sediment 
accumulation  rate  (w  =  3  cm  yr-!). 


NUC 

BUC  _ 

(mol  m~2  s“0 

NVC 

(mol  nr2  s-1) 

C:N  =  12.5 

C:N  =  54 

(mol  m~2  s“‘) 

Measured  SRR 
Model-calculated  R™- 
2KFc(iti>  (mol  nr2  sni')4 

Fmk  (mol  trr2  s~') 

XK/F 

1.5X10-7 

2.3X10”7 

2.7X10"* 

30 

4.8X10-7 

0.9X  10~7  4.5X10”7 

1.9X10“7  1.4X10-7 

9X1 0“8  (constant  for  all  stations) 

2.1  1.6 

9.6X10"7 

15.3X  10~7 
3.6X10-7 

4.0 

transport  of  the  solute  species  by  Spartina  root  systems  may 
not  be  quantitatively  significant,  compared  with  the  transport 
by  small  amount  of  polychaete  bioirrigation  in  the  vegetated 
portion  of  the  study  site. 

Although  phytoirrigation  appears  to  be  insignificant  as  a 
net  solute  transport  mechanism,  the  effect  of  vegetation  in 
the  salt  marsh  chemical  mass  transfer  is  important.  The  com¬ 
parison  of  IcG  profiles  between  NUC  and  NVC  (Fig.  2) 
makes  it  evident  that  the  vegetation  significantly  enhances 
the  net  OM  degradation  rate.  The  remineralization  of  OM 
occurs  predominantly  in  the  immediate  vicinity  of  WSI  at 
NUC  in  the  absence  of  Spartina  plants,  whereas  it  occurs 
throughout  the  sampling  depths  at  NVC.  Both  model-cal¬ 
culated  and  measured  SO reduction  rates  show  the  depth- 
dependent  decrease  in  SOj~  reduction  rates  at  NUC  and 
depth-independent,  quantitatively  significant  amount  of 
SOJ“  reduction  throughout  the  sampling  depth  at  NVC  (Fig. 
4).  It  can  be  speculated  that  Spartina  root  systems  enhance 
microbial  activities  by  modifying  the  habitats  through  the 
influx  of  02,  labile  OM,  and  nutrients,  as  well  as  through 
the  removal  of  toxic  metabolites  and  metals.  The  current 
modeling  has  limited  utility  in  phytoirrigation  study,  because 
it  does  not  explicitly  account  for  the  preferential  uptake  of 
certain  solute  species.  Such  an  account  would  require  a  bet¬ 
ter  understanding  of  the  chemical  mass  transfer  by  Spartina 
root  systems.  Further  investigations  of  the  concentrations, 
activity,  and  spatial  distribution  of  various  geochemical  spe¬ 
cies  and  microorganisms  in  relation  to  the  root  system  ge¬ 
ometry  are  warranted. 
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Fig.  6.  For  NUC,  an  elevated  C  :  N  ratio  (54)  yields  model-cal¬ 
culated  depth-concentration  profiles  that  agree  well  with  the  mea¬ 
sured  profiles. 


In  summary,  the  rates  of  OM  remineralization  and  TEA 
consumptions  were  determined  for  the  bioturbated  and  veg¬ 
etated  salt  marsh  sediments  of  the  Skidaway  SERF  using  the 
inverse  modeling.  The  models  were  constrained  by  the  eco¬ 
logical  data  as  well  as  measured  concentration  profiles  of 
major  redox  species,  and  the  model  results  were  verified  by 
the  measured  depth  profiles  of  SRR.  The  agreement  between 
measured  and  modeled  profiles  indicates  that  the  a  priori 
parameter  values,  including  those  associated  with  the  bio- 
turbation  and  irrigation,  are  reasonable.  The  calculated  rate 
profiles  not  only  reveal  the  difference  in  the  net  OM  deg¬ 
radation  rates  between  stations  but  also  verify  the  significant 
difference  in  the  relative  importance  of  Fe(III)- 
(oxy)hydroxides  and  SO^-  as  TEAs  between  stations  with 
varying  degrees  of  bioturbation.  In  heavily  bioturbated  BUC, 
Fe(III)-(oxy)hydroxides  are  recycled  >30  times  before  they 
are  permanently  buried  below  the  zone  of  bioturbation.  Re¬ 
cycling  at  NUC  and  NVC,  the  stations  that  host  few  biotur- 
bating  organisms,  is  limited  to  only  a  few  times.  The  model 
results  also  confirm  the  significance  of  vegetation  in  salt 
marsh  chemical  mass  transfer.  Although  the  phytoirrigation 
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Fig.  7.  For  NUC,  an  elevated  C  :  N  ratio  (54)  yields  model-cal¬ 
culated  SO2-  reduction  rate  profile  that  is  comparable  to  the  mea¬ 
sured  SRR. 
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of  major  solute  species  appears  to  be  insignificant,  Spartina 
root  systems  enhance  the  rate  of  OM  remineralization  in 
subsurface  sediments  well  away  from  WSI,  where  labile  OM 
and  microbial  activities  may  otherwise  be  limited. 
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